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Random Time-Scale Invariant Diffusion and Transport Coefficients 
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Single particle tracking of mRNA molecules and lipid granules in living cells shows that the time 
averaged mean squared displacement 5^ of individual particles remains a random variable while 
indicating that the particle motion is subdiffusive. We investigate this type of ergodicity breaking 
within the continuous time random walk model and show that 5^ differs from the corresponding 
ensemble average. In particular we derive the distribution for the fluctuations of the random variable 
5^. Similarly we quantify the response to a constant external field, revealing a generalization of the 
Einstein relation. Consequences for the interpretation of single molecule tracking data are discussed. 
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An ensemble of non interacting Brownian particles 
spreads according to Pick's law as a Gaussian packet. 
The ensemble averaged mean square displacement (MSD) 
is (x^(i)) = 2Dit where Di is the diffusion constant. By 
an Einstein relation Di is expressed in terms of statis- 
tical properties of the microscopic jumps according to 
Di = {Sx^) /2{t) where (r) is the average time between 
jumps and {Sx^) is the variance of the jump lengths. In- 
stead one can analyze the time series x{t) of the particle 
trajectory and determine the time averaged (TA) MSD 



52(A,i) = 



t-A 



(1) 



where A is called the lag time. For regular Brownian 
motion and long measurement time t 3> (t) we have 
S^ = 2DiA, i.e., an ergodic behavior such that the diffu- 
sion coefficient obtained from an individual trajectory is 
identical to the diffusion constant found from an ensem- 
ble of particles under identical physical conditions. 

From in vivo single particle tracking the diffusion of 
lipid granules in yeast cells [ Ij and of mRNA molecules 
in E.coli cells ^ two findings were made: (i) The TA 
MSD is subdiffusive, F ~ 2UaA" with a « 3/4. Usu- 
ally subdiffusion is defined by the behavior of an ensemble 
of particles {x'^{t)) = 2i:>„i"/r(l -I- a) and < a < 1. 
Such anomalous behavior is widespread |3|, [j, |5|, [y, |7|, 
including charge carrier transport in amorphous semi- 
conductors [J], models of gene regulation Jp, enzymatic 
binding in crowded cellular environments [9|] and anoma- 
lous dynamics of cell migration |10| to name but a few. 
(ii) The second striking observation [l|, l^] was that the 
TA diffusion coefficient D^ is a random variable different 
from the diffusion constant of the ensemble D^ , albeit the 
measurement time is long (see below) [ll|. Namely us- 
ing Eq. ([1]) to compute a TA MSD we get a result which 
varies from one single particle trajectory to another [l|, |2| 
(see Fig. 2]). This means that ergodicity is broken such 
that time and ensemble averages of the diffusion process 
are non identical. 

In this manuscript we investigate a widely applicable 
model for anomalous diffusion: the continuous time ran- 



dom walk (CTRW) 0, i, i, H, Q . In CTRW subdiffusion 
is scale invariant and (r) — > oo which naturally leads to 
ergodicity breaking [iJ, ll3( . We show that for the subdif- 
fusive CTRW the TA MSD ^ differs from the ensemble 
average, even in the limit of long averaging times. We ob- 
tain the distribution of TA MSDs that completely quanti- 
fies the magnitude of the new fluctuations. Then we treat 
the biased random walk showing that the TA response to 
an external driving field F also remains random. These 
new findings lead to a new type of fluctuation-dissipation 
relation for anomalous kinetics which depends both on 
the lag time A and the measurement time t. Finally we 
discuss the validity of our theory in experimental situ- 
ations and its generality in other models of anomalous 
diffusion. 

The uncoupled CTRW in one dimension is considered 
[i i, S, 0- The probability density function (PDF) 
of jump lengths is / (Sx) for which we assume that its 
variance (fe^) = J^^6x'^f{Sx)d{Sx) is finite. Waiting 
times between jump events are distributed with a com- 
mon PDF iP{t). So the particle waits in its initial loca- 
tion for a random waiting time, then makes a jump in 
space, and then the process is renewed. Our main inter- 
est is in the case where the average sojourn time is infinite 
(r) = oo namely the subdiffusive case with a power law 
PDF V(t) ~ Ar-(i+")/|r(-a)| and < a < 1. Phys- 
ical models which give specific values of a for different 
systems and models are given in Refs. [3|, y, [sj, y . 

We simulate CTRW trajectories with a = 3/4, for an 
unbiased random walk on a lattice /(x) — [d{x~l)+S{x+ 
l)]/2 and in Fig. [Jshow the TA MSD ^ of 10 individual 
trajectories with free boundary conditions. The most 
striking feature in the figure is that the curves are non- 
identical, and the TA MSDs remain a random variable 
even though a large number of jump events occur. In 
contrast if we choose a waiting time distribution with 
a > 1 the TAs will be identical to the ensemble average 
and non-random when the measurement time is long. 

To develop a theory for the observed behavior we first 
consider the average of Eq. H]) for unbiased CTRWs, 
namely for the case when the average jump length is zero. 
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FIG. 1: Simulations of the subdiffusive CTRW process with 
a — 3/4 and free boundary conditions show that the TA 
MSD is a random variable depending on individual trajecto- 
ries. The solid curve is the averaged behavior Eq. ((3}. The 
measurement time is t = 10* and ^(t) = aT~^^^°'' for t > 1. 



We consider first free boundary conditions, a widely ap- 
plicable case [3|, y, la, la, LZI , not necessarily relevant for 
the bounded motion in the cell (see below). We have 
x{t' -I- A) — x{t') — J2i=i'* ^^' ^^i where {6xi} are ran- 
dom jump lengths, and 'T.(t',t'+A) is the number of jumps 
in the interval {t',t' + A). For the unbiased CTRW the 
{Sxi}s are independent random variables with zero mean 
hence {[x(t' + A) - x(t')]^) = {Sx'^){n(^t,^t,+A))- The av- 
erage number of jumps in {t',t' + A) is (n^t'.t'+A)) = 
("(o,t'+A)) - ("(0,*'))- Using (uf^o.t')) - i'"/[Ar(l + a)] 



{S') 



{Sx^) fi 



it 



Ar(l + a) (l + a)(f-A) 

is obtained from Eq. ([1]). In the limit A ^ i we find 
,T7^ 2D„ A 



(6^) 



r(i + a)ii- 



(2) 



(3) 



where we used the generalized Einstein relation D^ = 
{Sx'^)/{2A) dil. For a 7^ 1 Eq. (13]) is very different 
from the behavior found for an ensemble, (x^(t)) — 
2Daf"/r(l -f a) indicating ergodicity breaking. Eq. ^ 
shows that if we know through measurement the ensem- 
ble averaged anomalous diffusion coefHcient Da we can 
determine the single particle trajectory averaged behav- 
ior. The result Eq. ([3]) can be explained by noting that 
the longer the process goes on the more likely we are to 
find long trapping times of the order of the measurement 
time (ageing). Hence (S'^) decreases when measurement 
time t is increased. Roughly speaking the diffusion con- 
stant depends on time D{t) ^ d{x'^)/dt ^ P^^ and Eq. 
([3]) is described by S'^ ~ 13 (t) A so a linear dependence on 
the lag time also seen in the simulations Fig. [I] is found. 



Distribution of 6^. For the CTRW under investigation 
we still have the usual scaling of cc^ ^ N with the number 
of jumps A'' in (0,i), however due to the broad distribu- 
tion of waiting times N ^ t" and so x^ ~ t". Similar 
scaling arguments can be used to analyze the distribution 
of 6'^. Assume no jump event occurs between time ti and 
t2 and that ^2 ^ ti+A. Then for ti < i' < i2 — A we have 
[x{t' + A) — x{t')]'^ — 0. Since for the scale free dynamics 
we have long sojourn times of the order of the measure- 
ment time without any jump event, [x{t'+A)—x{t')]'^ — 
for long renewal periods separated by shorter periods of 
activity. The most important point to realize is that for 
the process [a:(i'-|- A) —x{t')]'^ the distribution of sojourn 
times in state [x{t' -I- A) — x{t')]'^ = follows the same 
power law decay as the original process x{t) with a wait- 
ing time PDF tp{T) ^ r~(^+"^. This means that when N 
serves as the operational time wc have normal behavior 



(52 - CN/t 



(4) 



where C is a constant independent of N soon to be 
determined, and in the denominator we approximate 
t — A ~ f . Let PN{t) be the probability of making N 
jumps in the time interval (0,i) and Piq{u) its Laplace 
transform. From the convolution theorem Pm{u) = 
[1 — ■0(w)] exp[Af ln'i/;(u)]/u as well known [14]. Since we 
are interested in the long time behavior only the small u 
expansion tp{u) ~ 1 — Am" is relevant and we have 



Pn{u) - Au"-iexp(-7VAM"). 



Inverting to the time domain 

PN{t) 



aAl/aJ^l + l/a " I ^1/a^l/a 



t 



(5) 



(6) 



where la{t) is the one sided Levy stable PDF, whose 
Laplace pair is exp(— w") [15|,lla|. To find C we note that 
after averaging (P) = C{N)/t, using (N) - i"/yir(l+a) 
and Eq. ([3]) we have C = 2ADaA. By change of vari- 
ables we obtain the PDF of the dimensionless random 
variable ^ — P / {5"^) using Eqs. ([4]) and ([6|) 



hm (t>a (0 



ri/"(l + a) 



'a 



ri/"(l-|-a) 



_tl/Q 



(7) 



This is one of our main results since it describes the dis- 
tribution of a large class of time average observables, as 
we soon show. When a — > 1 we have an ergodic behav- 
ior \m\a^i(j)a{C) ~ '5(C ~ !)■ A measure of ergodicity 
breaking (EB) is the parameter 



EB = hm 



[(5^)')-{5^? 



(<52) 



2\2 



r(l + 2a) 



- 1 



(8) 



which is independent of the lag time A and Da ■ 

In Fig. [2] we show the behavior of the average of 5'^ 
and the fluctuations characterized by the EB parameter. 
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FIG. 2: (a) (52) versus A for a = 1/2 and f = 10*. Stars 
are simulations and the solid curve is theory Eq. ([Sj with- 
out fitting, (b) The EB parameter converges slowly to the 
asymptotic value EB = 0.5708 given by Eq. ((H}. Here A = 10 
(dots), A = 2500 (stars) and a = 1/2. 



FIG. 3: PDF of the scaled random variable ^ = S'^/{S^) for 
a = 1/2 and a = 3/4 with t = 10* and t = 10^ respectively. 
The full line is Eq. ^. Stars (A = 2500) and circles (A = 
4 X 10*) are simulations. 



showing excellent agreement between asymptotic theory 
and simulations though the convergence of the EB pa- 
rameter is typically slow. In Fig. [3] simulations of the 
PDF of (52/(^2^ foj. Q, ^ 2^/2 and a = 3/4 are_shown. We 
see that for a = 3/4 we have a peak close to 6'^/{S'^) = 1 
which indicates that we are closer to the ergodic phase 
(a — > 1) while for a — 1/2 the peak is on zero indicating 
stronger non-ergodic behavior as we decrease a. 

Biased CTRW and Generalized Einstein Relation. 
Now we assume that {Sx) ^ but constant, a case which 
leads to anomalous drift. We consider the TA 



S{A,t) 



t-A 



[x{t' + A) - x{t')]dt'/{t - A). (9) 



First we obtain the average which is done with an ap- 
proach similar to the unbiased case and using {x{t' + 
A)-x{t')) = {6x)[{t' + A)°' -t'°']/[Ar{l + a)] we find 



(S) 



{5x) 



A 



AT{l + a)t^- 



(10) 



for t > A. Then we can show that the PDF of ^ = 5/ {5) 
is given by Eq. (0 thus fluctuations of the TA MSD 
of the unbiased random walk and the fluctuations of the 
biased mean response have identical distributions. 

As is well known according to the generalized Ein- 
stein relation the transport of an ensemble of particles 
is related to the free diffusion of the same particles by 
{x{t))F = F{x^{t))/2kbT where {x'^{t)) is the ensemble 
average MSD in the absence of a force field and {x{t))p 
is the mean drift when a constant force F is applied to 
the system [3|, |l7|, UM, uM ■ This relation can be used to 
prove that on a microscopic scale {Sx)f = F{6'^x)/{2kbT) 



which can be obtained from thermal detailed balance con- 
ditions JJ}. Using this relation and Eqs. ^ and (fTUl) 



{d)F = F{S^)/[2hT] 



(11) 



which is valid under usual linear response assumptions. 
This Einstein relation for the TAs while clearly related 
to the Einstein relation for the ensemble average, is valid 
for any lag time A and measurement time t. In this sense 
it differs from the usual Einstein relation. As mentioned, 
the relation between transport and diffusion runs deeper, 
at least within the CTRW model, since we showed that 
the fluctuations are identical as long as the external field 
does not modify iP{t) (a reasonable assumption for weak 
fields TtI) and described by Eq. ^. 

Relation with experiments. Our results can be tested 
in single particle experiments for example for a bead 
anomalously diffusing in an actin network which exhibits 
a CTRW type of dynamics [20|. However in the exper- 
iments in the cell [l|, l2] the particle motion is bounded 
by the cell walls. Indeed the particles may interact with 
the cell wall many times whenever 2Dat°'/r{l + a) > L^, 
where L is the system length. Finiteness of the system 
implies that at long times the ensemble averaged MSD 
will not increase with time but rather saturate. 

We have simulated the effect of a boundary by consid- 
ering an unbiased CTRW on a lattice with system size 
L = 62 with lattice spacing equal unity and a = 3/4. 
The simulations shown in Fig. [4] look similar to experi- 
ment and we have S"^ ~ A'^ with (3 ~ 3/4 at least within 
a reasonable time window. The exponent /3 < 1 depends 
on the system size, on a and on the time window un- 
der investigation. Still our numerical results show that 
CTRW theory is compatible with available experiment. 
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nents is described by Eq. ([7]); similar behavior is found 
for super-diffusive Levy walks. For random walks in ran- 
dom environments (e.g., random trap and comb models) 
we expect similar behavior due to the deep connections 
between these models and CTRW theory (in these models 
(r) ^ cx) as in the subdifFusive CTRW however the disor- 
der is quenched not annealed) . Thus one of the most ba- 
sic paradigms of transport and diffusion theory, namely 
that information obtained from single particle tracking 
is contained already in the ensemble measurement is not 
valid for anomalous diffusion. This has ramifications for 
vast classes of processes. 

This work was supported by the Israel Science Foun- 
dation and the Kort fellowship. We thank Ido Golding 
and Lene Oddershede for discussions. 



FIG. 4: Time average 5^ vs. A for unbiased CTRW on a 
lattice of size L = 62, mimicking a particle bounded in a finite 
domain, as found in the cell. Compared with the unbounded 
case in Fig. [T]the diffusion is slower. The trajectory averaged 
MSD follows 52 ~ A'' with /3 = 3/4 similar to what is found 
in [3, [2]. The measurement time was t — 10*, a = 3/4. 



A direct test of our theory would be to change the ex- 
perimental time t (not only A as done so far) and see if 
the TA diffusion slows down with increasing t. 

Notice that for free boundary conditions we get 5'^ ex 
A, the existence of the boundaries thus causes the dif- 
fusion to appear slower (i.e. /? < 1) which is intuitively 
expected. For free boundary conditions what appears as 
normal diffusion in a single particle measurement there- 
fore may actually be a hidden subdiffusive process. For 
both free and reflecting boundary conditions as we in- 
crease the measurement time diffusion is slowed down 
when the TA procedure is made. In all such experiments 
it is thus imperative to analyze the TA MSD also as func- 
tion of the measurement time t. Additional clues about 
the nature of the diffusion are the potential scatter of the 
diffusivity as well as the shape of the trajectories. 

In the cell, the measurement time t is limited by the 
life time of the cell. This is important from the point 
of view of theory which usually assumes an infinite mea- 
surement time. Indeed in a finite volume one can expect 
from a thermodynamical argument demanding station- 
arity that if ipir) decays like a power law it does so only 
within a finite time interval and then a cutoff will appear. 
However the finite life time of the cell implies that the 
usual long time limit essential for ergodicity may not be 
reached and ergodicity breaking is found: Ergodicity of 
diffusion processes is not fulfilled in a living cell. 

More generally we expect TA diffusion and transport 
coefficients to remain random in other models of anoma- 
lous diffusion. We have recently shown 2l| that for in- 
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